Genome-Wide Association Study (GWAS) and genome prediction of seedling salt tolerance in bread wheat (Triticum aestivum L.)

Background Salinity tolerance in wheat is imperative for improving crop genetic capacity in response to the expanding phenomenon of soil salinization. However, little is known about the genetic foundation underlying salinity tolerance at the seedling growth stage of wheat. Herein, a GWAS analysis was carried out by the random-SNP-effect mixed linear model (mrMLM) multi-locus model to uncover candidate genes responsible for salt tolerance at the seedling stage in 298 Iranian bread wheat accessions, including 208 landraces and 90 cultivars. Results A total of 29 functional marker-trait associations (MTAs) were detected under salinity, 100 mM NaCl (sodium chloride). Of these, seven single nucleotide polymorphisms (SNPs) including rs54146, rs257, rs37983, rs18682, rs55629, rs15183, and rs63185 with R2 ≥ 10% were found to be linked with relative water content, root fresh weight, root dry weight, root volume, shoot high, proline, and shoot potassium (K+), respectively. Further, a total of 27 candidate genes were functionally annotated to be involved in response to the saline environment. Most of these genes have key roles in photosynthesis, response to abscisic acid, cell redox homeostasis, sucrose and carbohydrate metabolism, ubiquitination, transmembrane transport, chromatin silencing, and some genes harbored unknown functions that all together may respond to salinity as a complex network. For genomic prediction (GP), the genomic best linear unbiased prediction (GBLUP) model reflected genetic effects better than both bayesian ridge regression (BRR) and ridge regression-best linear unbiased prediction (RRBLUP), suggesting GBLUP as a favorable tool for wheat genomic selection. Conclusion The SNPs and candidate genes identified in the current work can be used potentially for developing salt-tolerant varieties at the seedling growth stage by marker-assisted selection. Supplementary Information The online version contains supplementary material available at 10.1186/s12870-022-03936-8.

worldwide is salinity/salt stress, which damages numerous physiological, biochemical, and molecular processes [5,6]. Salinity is one of the important abiotic stresses that can seriously disrupt wheat production [7]. Generally speaking, when neutral soluble salts (chlorine, calcium, sodium, etc.) excessively accumulate in the rhizosphere, they can disrupt nutrient uptake [8]. Therefore, excess salts in the soil can lead to nutrient imbalance and ionic toxicity/deficiency, which negatively affect wheat yield [8][9][10]. Previous studies have demonstrated yield loss of up to 50% in wheat when exposed to a high salt concentration [11]. Thus, there is a demand to uncover salinity-responsive genes and use them to develop new salt-tolerant varieties [12].
Salt tolerance is a complex trait that includes a variety of genes, regulation networks, signal transductions, and metabolic pathways [13][14][15][16]. On the other, wheat response to saline environments depends on the duration and intensity of the stress and differs between genotypes as well as growth stages [17,18]. For these reasons, assessing a genetic panel for salt tolerance at the seedling growth stage is a difficult task for wheat breeders. To make further progress in the development of salinity-tolerant wheat varieties it is crucial to get a better understanding of the molecular basis of salinity tolerance-related traits by using genetic tools, like quantitative trait loci (QTL) mapping [19].
QTL mapping has been used for detecting genes/ genomic regions linked to salt tolerance traits, such as bio-physiological (e.g., Na + /K + ratio) and agronomical traits (e.g., grain yield) in the salt-stressed wheat fields [7,15,19]. Importantly, these endeavors have relied on mapping populations of small size and a low number of SSRs markers, reflecting a limited resolution of QTLs, which cannot be reliably adopted in the marker-assisted selection. In contrast, genome-wide association study (GWAS) provides an alternative to QTL mapping for identifying genes linked to the phenotype of interest [20]. Association mapping can be performed by single-locus (GLM and MLM) or multi-locus (mrMLM) models [21]. GLM and MLM models adopt a genome scan by testing SNP markers at a time and need strict multiple test correction (e.g., Bonferroni) for managing false positives. However, this process is often too conservative and may lead to the loss of statistical power, failing to detect true associations that may be important. Moreover, singlelocus models cannot simultaneously estimate all marker effects, and thereby cannot present a proper model for complex traits, which are controlled by the cumulative effect of several genes. To overcome these challenges, multi-locus approaches have started to be widely adopted as an alternative approach for dissecting the molecular basis of quantitative traits in plants and crops .
Previous studies presented experimental evidence regarding the QTLs/ candidate genes related to the salt tolerance at the seedling stage (i.e., seedling salt tolerance) in various plants/crops. In a research attempt, Luo et al. [30] elucidated the genetic basis of seedling salt tolerance by 557,894 polymorphic SNPs on 348 maize inbred lines. They identified 13 candidate genes associated with seedling salt tolerance by GWAS, among which, ZmPMP3 and ZmCLCg were confirmed as genes involved in seedling salt tolerance. Interestingly, ZmCLCg was found as a chloride transport in maize. By using 18,430 polymorphic SNPs on 149 cotton genotypes, Zheng et al. [7] found six seedling salt tolerance genes, including Gh_D08G1309, Gh_D08G1308, Gh_A01G0908, Gh_A01G0906, Gh_D01G0945, and Gh_D01G0943, which were found to be responsible for cell amplification, auxin response, N-glycosylation, transmembrane transport, osmotic pressure balance, sucrose synthesis, and intracellular transport, respectively. Thabet et al. [28] evaluated 121 barley accessions for seedling salt tolerance by using 9 K SNPs and revealed around 1500 candidate genes, which encode potassium channels mapped on Ch.1H. The squamosa promoter-binding-like protein 6 at Ch.5H was detected to be linked with seedling salt tolerance. Screening a total of 203 rice accessions led to uncovering of 26 QTLs for seedling salt tolerance. Candidate genes for promising QTLs included glycosyl hydrolase, sucrose transporter, leucine zipper TF, ammonium transporter, and MYB TF [48].
As auxiliary tools for GWAS, genomic prediction boosts the speed and effectiveness of breeding by decreasing the time required for breeding cycles and by increasing selection accuracy [23]. The marker set, genomic selection method, population structure, and trait genetic architecture are the main factors that impact genomic accuracy. Several projects have demonstrated moderate to high genomic accuracy for complex traits in barley [24], maize [25], oat [26], rice [27], and wheat [23]. However, genomic prediction of the salt tolerance at the seedling stage has not been reported in wheat.
To the best of our knowledge, little is known about genomic regions associated with salt tolerance at the seedling stage in wheat. Therefore, we uncovered putative candidate genes and evaluated the genomic prediction accuracy of salt tolerance at the seedling stage using three methods for building a genomic selection model, namely GBLUP, RRBLUP, and BRR.

Traits phenotyping
The phenotypic evaluation showed that most seedlingrelated traits have lower performance under salinity than normal conditions, highlighting salt stress limits

Linkage disequilibrium (LD)
In the panel of cultivars, LD calculation using 46,203 SNPs led to the detecting of 1,830,925 markers pairs (MPs), of which 60% of them displayed significant linkage. LD between marker pairs was recorded across the 21 chromosomes ranging from 0.14 (Ch.6D) to 0.37 (Ch.4A). The highest number of MPs were discovered in the B genome (949,425, 51.85%), followed by the A genome (675,325, 37%) and D genome (206,175, 11.26%) ( Table 4).
Implementing a similar test on wheat landraces led to uncovering 1,828,675 MPs with a mean r 2 of 0.18, which is lower than that in wheat cultivars. Of course, a bigger part of marker pairs was found significant (836,400, 45.74%) in landraces. LD was strongest between marker pairs in Ch.4A (0.32), followed by Ch.2A (0.25) ( Table 4).

Population kinship and structure matrix
Based on the ∆K formula, the optimum number of subpopulations (K) in the association panel was estimated at K = 3 (Fig. 2S). From the PCA, first two PCs explained 17.0 and 6.4% of the genotypic variance, respectively (Fig. 1). Clear subpopulations were observed from the first two PCs, which indicated three subpopulations with admix accessions falling between clusters. As the panel of wheat cultivars and landrace have subpopulations, the PCA and kinship matrix were performed as variancecovariance. The cluster analysis based on the kinship matrix exhibited that the SBP-I subpopulation harbors 110 accessions (105 landraces and 5 cultivars), the SBP-II harbors 38 accessions (28 landraces and 10 cultivars), and the SBP-III harbors 144 accessions (69 landraces and 75 cultivars) (Fig. 2). A neighbor-joining tree of all accessions also clearly exhibited the clustering into three subgroups (Fig. 3).

MTAs for seedling-related traits
Using mrMLM model, 817 and 1006 significant MTAs were identified under normal and stress conditions, respectively, for morphological, physiological, and biochemical traits at -log 10 (P) > 3 (Fig. 4). Among these, 40 and 29 highly significant, functional MTAs were regarded as "reliable" MTAs under normal and stress conditions, respectively. The reliable MTAs were selected based on the fact that they passed a high significance threshold and also have a cellular function. From the reliable MTAs, we selected "major" MTAs, which explained ≥10% of the phenotypic diversity for the traits. A total of 15 and 8 major MTAs were detected for control and salt stress, respectively (Tables 5 and 6). QQ and Manhattan plots of top SNPs for the traits of interest are presented in Fig. 5.

Putative candidate genes for salt tolerance
The analysis of gene ontology on 29 reliable MTAs indicated that the candidate genes harboring these SNPs encode proteins involved in several biological processes, including photosynthesis, response to abscisic acid, cell redox homeostasis, sucrose and carbohydrate metabolism, ubiquitination, transmembrane transport, and chromatin silencing under salt stress. From the homologs in rice (Tables 7 and 8), 25 putative candidate genes were detected for response to salt stress.

Genomic prediction (GP)
Under stress, the highest genomic prediction accuracy was achieved for RWC, ELI, chlorophyll, carotenoid, protein, and CAT traits by the GBLUP method. By the RR-BLUP method, the highest prediction accuracy was observed for GPX, root volume, and K + content traits. The BRR method showed the highest prediction accuracy for SPAD and proline traits (Fig. 6). Overall, the GBLUP model exhibited better performance than BRR and RR-BLUP, suggesting that GBLUP is the preferable tool to use for genomic selection in the wheat panel.

Discussion
Breeding for salt tolerance in wheat is a challenging task due to the polygenic nature of this trait and the polyploid nature of the wheat genome. This task is further complicated by the fact that various mechanisms are adopted for salinity tolerance at the seedling and adult growth stages [24]. To the best of our knowledge, little is known about genomic regions associated with salt tolerance at the seedling stage in wheat.
With such a situation in mind, we developed a GWAS panel consisting of 298 Iranian bread wheat accessions and used this panel to identify candidate genes involved in controlling salinity tolerance at the seedling stage.

The impact of salinity on wheat seedling traits
In-depth phenotyping is a key part of a GWAS procedure [29]. Herein, a total of 25 seedling-linked traits were evaluated that have been previously employed for QTL mapping of salinity tolerance at the seedling stage in cotton, rice, and maize [7,9,10]. Similar to our observations, previous reports have also shown that salinity negatively affects seedling-related traits [29][30][31][32]. In a conclusion, salt stress remarkably limits wheat seedling growth, as previously reported by Liang et al. [9]. From our findings, a negative correlation was found between Na + levels and root volume, showing the detrimental effect of sodium ions on the root system. The inherent capability of accessions to maintain low Na + levels is thus one of the critical parameters inducing salt tolerance. Other mechanisms for salt tolerance include tissue tolerance and Na + compartmentalization which may be also involved in salinity tolerance at the seedling stage in wheat accessions [33].

Population structure of the wheat panel
Structure analysis disclosed three subpopulations among 298 Iranian bread wheat accessions. The results from   the PCA also support this observation. Interestingly, the clustered pattern of wheat accessions was not consistent with their geographical distribution or origins (Table S1,  Table S2, and Fig. 3). This can be likely attributed to the migration of farmers to different regions and germplasm exchange across institutes and researchers across the world [32].

Linkage disequilibrium in wheat sub-genomes
In line with previous reports, most markers were located in the B and A genomes [34], and the same trend was recorded for MPs in LD. The higher variation observed in the A and B genomes is likely a consequence of two factors [35], the older evolutionary history of these genomes and gene flow from the species T. turgidum (but not Ae. tauschii) to common wheat. From our observations, LD and marker distance across the A and B genomes were much lower than in the D genome. The fact that cultivars exhibit higher LD in contrast to landraces is likely a result of selection events during crop breeding [23]. In addition to selective breeding, other factors such as recombination, population relatedness, genetic drift, mutation, and mating systems affect linkage disequilibrium in wheat and other plants [36].

Candidate genes for salt tolerance at the seedling stage
To date, many genes and QTLs connected with salinity tolerance at the seedling stage have been reported by association and linkage mapping in various crops and plants. However, little is known about the link between genomic regions associated with seedling salt tolerance with corresponding mechanisms in bread wheat. We successfully identified 27 putative candidate genes for salinity response that encode proteins/ enzymes involved in antiporter, electron transfer, kinase, hydrolase, endoribonuclease, ATPase, glutamate receptor, metalloaminopeptidase, glycosyltransferase, oxidoreductase, acyltransferase, calcium ion binding, ubiquitin transferase, sucrose synthase, etc. From mapping wheat SNPs on the rice genome, 25 putative candidate genes, including OsPAP1d, OsPAP1c, OsIDI4, OsGPCR, OsENODL6, OsGELP83, OsWD40, OsRFPH2, and OsRLCK202 were shown to be responsive to salinity. We must remind that the genomic regions associated with seedling salt tolerance, it is a problematic comparison across various studies because of the difference in the mapping population and marker platforms, as well as the absence of a consensus map for comparing genomic locations.

Candidate genes for root/shoot height and weight
Root and shoot height and weight are key traits that specify plant architecture and affect grain yield in salt environments. The genetic basis of these traits is complex, and controlled by many genes and the environment [32]. To date, several genes have been found to be responsible for controlling root/shoot height and weight at the seedling stage of various plants [10,[28][29][30][31][32]. In this study, the markers rs53540, rs35884, rs257, rs37983, rs18682, rs55629, and rs44076 were linked to shoot fresh weight, shoot dry weight, root fresh weight, root dry weight, root volume, root length, and shoot height traits, respectively, allowing the identification of reliable salt-responsive genes. Among these, TraesCS1D02G156100, TraesCS3B02G182700, TraesCS7B02G339500, TraesCS3B02G227800, TraesC-S4A02G415700, and TraesCS1B02G480700 explained a large fraction of the phenotypic variance (≥ 10%) and classified as "major" candidate genes. Which can be targeted in future research. From mapping, the wheat SNPs on the rice genome, the root volume-connected SNP on the rice Ch.9 led to the detecting the IDI4 gene of 1-aminocyclopropane-1-carboxylate synthases family, which have a critical function in response to hypoxic stress in crops [37].

Candidate genes for RWC and proline content
Two major candidate genes TraesCS1D02G156100 and TraesCS4A02G415700 were identified that control RWC and proline and are located on Ch.1D and Ch.4A, respectively. From mapping the wheat SNPs on the rice genome, one proline-related SNP on the rice Ch.7 led to discover of a member of the WD40 protein family, WD40-145, which response to salt stress likely through interaction with MADS-box, MYB, and bHLH TFs [38]. Interestingly, the SPAD-connected SNP on the rice Ch.11 revealed a 2,3-oxidosqualene cyclase (OSC7), which constructs the skeleton of cyclic triterpenoids [39]. Terpenoids produced by oxidosqualene cyclases, such as α-or β-amyrin, play an essential role to cope plant roots with salinity [40].

Candidate genes for CAT and GPX activities
In the salt-stressed seedlings, the rs10254 and rs61179 markers were detected to be associated with CAT and GPX activities, highlighting the effect of the reliable responsive genes TraesCS3B02G556500 and TraesCS1B02G048300, respectively. From mapping the wheat SNPs on the rice genome, the homolog genes Os05g0121900 and Os07g0105600 were uncovered for affecting CAT and GPX activities on the rice Ch.5 and Ch.7, respectively. The former codes a phosphate/phosphoenolpyruvate translocator (PPT) protein-like, which is responsible for the development of phenylpropanoid metabolism-derived signal molecules triggering leaf intervene regions [41], and the latter codes a photosystem II oxygen-evolving complex protein, which is involved in transferring electrons within the cyclic electron transport pathway of photosynthesis.

Candidate genes for pigment contents
Salt stress can inhibit PSII activity and destroy chlorophyll molecules, ultimately influencing a plant's ability to photosynthesize [38]. To date, several QTLs for chlorophyll content has been identified during early growth stages under salinity. In our experiment, markers rs34693, rs18445, rs34693, and rs59624 were associated with to chlorophyll a, chlorophyll b, total chlorophyll, and carotenoid traits, highlighting the reliable responsive genes TraesCS7B02G289500, TraesCS6A02G347900, TraesCS7B02G289500, and TraesCS6B02G343300, respectively. Interestingly, the homolog gene CYP97A4 was earlier identified as it influenced chlorophyll b content. Similarly, Chaurasia et al. [33] identified a gene encoding cytochrome 450, CYP709B2, which was involved in regulating leaf chlorophyll levels. CYPs are known to play a key role in response to salt stress by hormone signaling and/or through accelerating ROSs scavenging. Kushiro et al. [25] also uncovered an Arabidopsis CYP gene, CYP709B3, which is responsible for ABA signaling and salt response. Overall, our observation suggests that the CYP gene identified from the chlorophyll-related SNP may have a vital function in specifying wheat response to saline soils. Le et al. [43] found two SNPs for chlorophyll content located in the genes OsRLCK253 (Ch. 8) and OsCYL4 (Ch. 9) in saltstressed rice. The first gene encodes a receptor-like kinase, which is known to be involved in salinity tolerance, while the second code a cyclase-containing protein, which negatively regulates stress tolerance linked to ROS levels. Le et al. [43] also detected several genes associated with chlorophyll b content, including OsNUC1 (Nucleolin-like protein), OsHox33 (HDZIP III TF), OsARF25 (Auxin response factor), OsWAK128 (OsWAK receptor-like kinase), OsCHX15 (ATCHX protein), and OsZFP213 (C2H2 TF). Moreover, we discovered one MTA for total chlorophyll content that was linked to OsENODL6 homolog, which encodes an early nodulinlike protein in rice (located on Ch.2). Early nodulin-like proteins have been shown to display ≥3-fold changes in salt-stressed Cajanus cajan plants, thus, Awana et al. [42] suggested their involvement in the salt response.

Candidate genes for pigment contents
From earlier studies, genotypes tolerant to saline environments can decrease osmotic stress, absorb more K + , and prevent Na + accumulation in order to maintain a low Na + /K + ratio [33]. Thus, Na + and K + -related genes were explored in our experiment to figure out K + and Na + -dependent wheat responses to salt stress at the seeding stage. In a high salt environment, Na + toxicity and osmotic imbalance are two limiting factors for crop growth [12,33]; so researchers have linked Na + exclusion capability to grain yield under salinity stress [11]. Therefore, genes related to low Na + content are key candidates for improving salt tolerance in wheat. Earlier studies have detected genomic regions associated with Na + exclusion on Ch. 1A, 2A, 2B, 5B, and 6B in salt-stressed wheat [16]. Interestingly, we uncovered TraesCS1B02G472200 and TraesCS4B02G330600 as genes associated with Na + accumulation in the shoot and root, respectively, suggesting these genes may play significant roles in sodium homeostasis at the wheat seedling stage. Chaurasia et al. [33] found three major QTNs for Na + content in wheat (Q.Na-6DL, Q.Na-6AL, and Q.Na-2AS), among them, Q.Na-6DL had a remarkable contribution to Na + accumulation. From mapping the wheat SNPs on the rice genome, the root Na + content-related SNP on the rice Ch.4 led to the detecting of a member of RFPH protein family, OsRFPH2-14, which operates as RING-H2 Finger E3 ubiquitin ligase. Similarly, Liu et al. [45] reported that the OsRFPH2-10 gene reduces the level of P2 protein and incorporates antiviral defense at the early infection stage.
In addition to Na + , K + homeostasis is important for crop tolerance to salinity, since this ion is responsible for many key physiological processes like stomata movement, protein synthesis, respiration, photosynthesis, and growth metabolic functions [46]. In fact, higher K + content may enable wheat to tolerate salt stress by developing a root system. We successfully identified TraesC-SU02G075800 and TraesCS5A02G109600 as genes linked with K + concentration in the shoot and root, respectively, suggesting these genes are important for K + homeostasis at the wheat seedling stage. From the mapping of wheat SNPs on the rice genome, the root K + content-related SNP on the rice Ch.6 revealed the receptor-like cytoplasmic kinase 202, OsRLCK202. Differential expression patterns of OsRLCKs at various development stages and stress suggest its involvement in diverse functions. Lin et al. [47] found a genomic region on Ch.1 associated with shoot K + content (OsHKT1) that explained 40% of the phenotypic variation. Map-based cloning showed that this gene encodes a Na + transporter, HTK1, which is responsible for K + and Na + homeostasis.
The K + /Na + ratio is a well-known index that reflects a whole-plant response to salt stress. Generally speaking, salinity-tolerant accessions hold a low ratio of Na + / K + in aerial parts [48]. Genomic regions related to this trait have been detected in different plants and crops and   attempts are currently being made to use them in the development of high-yield cultivars tolerant to saline soils [16]. Earlier studies have reported the genomic regions on 2AL, 4AS, and 7DL associated with Na + /K + ratio in saline fields [33,47]. We successfully identified the genes TraesCSU02G082000 and TraesCS6D02G403800 for K + /Na + ratio in shoot and root, respectively, indicating potential targets for salt tolerance breeding. Chaurasia et al. [33] reported a novel QTN (Q.NaK-1BS) for K + / Na + ratio on 1BS in wheat that explain 4-38% of the phenotypic variation. Annotation of this locus demonstrated that Q.NaK-1BS is located inside the Rab-like-GTPase gene, which plays a vital function in salt tolerance by regulating Na + transportation [49]. Batayeva et al. [48] found one genomic region associated with the Na + /K + ratio on rice Ch.3 that harbored a sucrose transporter gene. Finally, Li et al. [50] discovered one novel QTL (qSNK3-1) located on rice Ch.3 that explains 14% of phenotypic variation. This QTL coincided with OsIRO3 gene, which encodes a bHLH-type TF and acts as an inhibitor of Fe-deficiency response in rice.

Genomic selection in wheat panel
The GP accuracy depends on the genomic selection method, level of LD, genetic diversity in the studied population, and genetic architecture of the studied trait [23]. In this study, we observed that the GBLUP method had better performance than the RR-BLUP and BRR methods, suggesting that GBLUP is a powerful tool for implementing genomic selection in wheat. Previous studies have suggested that high prediction accuracy can be achieved by GBLUP if markers are closely linked to the trait of interest. RR-BLUP works well for traits where the genetic architecture consists of numerous loci with small effects while the BRR approach is similar to RR-BLUP, except marker effect shrinkage depends on population size in BRR [23]. The better performance of GBLUP in our study could depend on the fact that SNPs in this study were closely associated with salt tolerance traits at the seedling stage in wheat.

Conclusion
Our work provides new insights into the molecular mechanisms underlying salt tolerance traits at the seedling stage in wheat. Putative candidate genes controlling these traits, i.e. K + /Na + ratio, can be targeted for developing salt-tolerant wheat cultivars at the seeding stage using marker-assisted selection. Moreover, genomic selection by using our putative genetic markers along with GBLUP-based genomic prediction will help to achieve the above-mentioned goal. Identification of varieties with high salt tolerance at the seedling stage, as well as knowledge of the associated SNPs and haplotype, could be useful for wheat production and for improvement of direct-seeding varieties.

Plant material
A total of 298 Iranian bread wheat genotypes were evaluated in this study. The wheat panel contained 90 cultivars released during 1942-2014 and 208 landraces gathered from the Persian plateau during 1931-1968. All the materials were provided by the Seed and Plant Improvement Institute and the Tehran University, Karaj, Iran. More details on these bread wheat accessions can be found in Tables S1 and S2.

Experimental design and phenotyping
The wheat cultivars and landraces were assessed for salt tolerance at the seedling stage using two salinity levels: 0 (control) and 100 (stress) mM NaCl (the selection of 100 mM NaCl stress was based on previous studies and the tolerance threshold of wheat to salinity). The study was carried out in a factorial experiment-completely randomized design (CRD) with two repeats and two factors: the first factor accounting for 298 Iranian bread wheat accessions and the second factor for two salinity concentrations. For each treatment, eight healthy and surface-sterilized seeds from each accession were planted in plastic pots (2 kg, 14 cm diameter, and 14 cm height). The soil composition of each pot was made up of a 3:2:1 ratio of decomposed litter, soil, and sand, respectively. The average temperature in the greenhouse was set to 25 °C during the day and 20 °C during the night, with a 6 h light/8 h dark photoperiod and 60% relative humidity. A thinning step was carried out at the two-leaf stage and four seedlings remained in each pot. Salt stress was gradually applied 15 days after germination by adding NaCl (25 mM) every other day together with irrigation water to reach the final concentration of NaCl, i.e., 100 mM. Crops Fig. 6 The impact of genomic selection (GS) methods on genomic prediction (GP) accuracy for 25 various traits in Iranian wheat landraces and cultivars under normal and salinity conditions. The prediction accuracy for RR-BLUP, GBLUP, and BRR-based genomic selection (GS) is presented with green, red, and blue colors, respectively. The middle point of boxplots indicates a mean of GP accuracies for the trait of interest. Electrolyte leakage (ELI); SPAD; shoot fresh weight (SFW); shoot dry weight (SDW); relative water content (RWC); root fresh weight (RFW); root dry weight (RDW); root volume (RV); shoot height (SH); root height (RH); chlorophyll a (Chl a); chlorophyll b (Chl b); total chlorophyll (total Chl); carotenoid (Car); protein; proline; catalase (CAT); guaiacol peroxidase (GPX); malondialdehyde (MDA); Shoot Na (Na-s); Root Na (Na-r); Shoot K (K-s); Root K (K-r); Shoot K/Na (K/Na-s); root K/Na (K/Na-r) were harvested three weeks after stress to measure the following morpho-physiological characteristics with two repeats: root volume (RV), root length (RL), shoot height (SH), root dry weight (RDW), shoot dry weight (SDW), root fresh weight (RFW), shoot fresh weight (SFW), malondialdehyde (MDA), electrolyte leakage (EL), relative water content (RWC), proline (P), soluble protein (PC), catalase (CAT), guaiacol peroxidase (GPX), photosynthetic pigments, SPAD, Na + content, K + content, and K + /Na + ratio.

Physiological trait measurements Electrolyte leakage (EL)
Identical circular pieces were prepared from fully-developed leaves and placed separately in plastic-capped tubes containing distilled water for 24 h at room temperature after which the solution's electrical conductivity (EC 1 ) was measured. The tubes were put in a Ben Marie apparatus at 95 °C for 90 min, and after cooling to 25 °C, electrical conductivity (EC 2 ) was measured. The EL% was calculated as (EC 1 / EC 2 ) × 100.

Leaf greenness
This trait was evaluated by using a SPAD-502 plus chlorophyll meter. Greenness levels were recorded based on the mean of three sections from the youngest fully-developed leaves.

Relative water content (RWC)
The highest leaves were harvested and their fresh weights (FW) were measured immediately. To determine the turgid weights (TW), the leaves were put down in distilled water overnight at low light intensity (to limit weight loss due to respiratory activity) and then weighted again. Eventually, leaves were placed at 70 °C for 48 h and their dry weights (DW) were recorded. Relative water content (%RWC) was estimated as: [(FW-DW)/(TW-DW)] × 100.

Proline content
Proline level was measured using the method developed by Bates et al. [51]. Briefly, 0.5 g of the fresh leaf was mixed with 10 mL of 3% sulfosalicylic acid and completely homogenized in a mortar. To remove excess materials from the solution, the tubes were centrifuged for 15 min at 15,000 rpm, 4 °C. The solution (2 ml) was mixed with 2 mL of ninhydrin and 2 mL of acetic acid. The tubes were kept in a hot water bath for 1 h and then cooled down in an ice bath for 1 h. Tubes containing 4 mL of toluene were vortexed for 20s and the proline content of the supernatant was estimated by a spectrophotometer at 520 nm.

Total protein
Leaf protein content was estimated based on Bradford [52]. Briefly, 500 mg of fresh tissue was homogenized in 5 mL of potassium phosphate buffer (10 mM, pH 7) with 5% (w/v) PVP, followed by centrifuging for 25 min at 15,000 rpm, 4 °C. Bradford reagent (990 μL) was mixed with 25 μL of supernatant and absorbance was read at 595 nm.

Malondialdehyde (MDA)
To detect MDA levels, as an output of lipid peroxidation, the plant extract was prepared using 1.0 g of tissue as explained by Cakmak and Horst [53]. After recording absorbance at 600 and 532 nm, the 155 mM − 1 cm − 1 extinction coefficient was used in the following formula to estimate the MDA level: nM MDA = A 532 -A 600 /1.55*10 5 .

Antioxidant enzyme activities
To prepare the enzymatic extract, 0.1 g of fresh tissue was crushed in liquid nitrogen, followed by adding 1 mL of sodium phosphate buffer (50 mM, pH = 7). The homogenate was centrifuged for 20 min at 10,000 rpm and 5 °C after which the CAT and GPX activities were measured from the resulting supernatant [53]. The enzyme activities were expressed as changes in absorption/min/g of fresh weight.

Photosynthetic pigments
Carotenoid and chlorophyll (a, b, and total) levels were measured based on the procedure described in Arnon [54]. Light absorption was read at 645 and 663 nm by a spectrophotometer and the chlorophyll levels were determined as follows: Where A is the optical absorption of samples, V is the ultimate acetone volume, and W is the leaf fresh weight.
The total carotenoid was calculated as follows: A 1% 1cm × 100 × W K + /Na + ratio, Na + content, and K + content Three leaves of individual accessions were gathered and dried for 3 days at 55 °C and 0.5 g of dried leaves were cut into pieces and put in a digestion tube (100 ml). A total volume of 10 mL of HClO 4 and HNO 3 (at a 1:3 ratio) was added to the tubes. The tube was then put in a digestion block for heating for 2 days. After cooling the transparent extract, the flasks were calibrated to a final volume of 25 mL by adding distilled water. By using a Flame Photometer, the K + and Na + contents were estimated from the filtered solution [55].

Phenotypic data analysis
The variance analysis (ANOVA) of data collected in the normal and salinity environments was implemented by SAS 9.4 (SAS Institute, USA). The analysis was followed by calculating Pearson's correlation coefficient to disclose significant relationships (P < 0.01) between traits. The descriptive statistics of phenotypic datasets were calculated by SPSS Statistics 21.0 (IBM Inc., USA).

Genotyping and SNP imputation
The genomic DNA was extracted from wheat seedlings by the CTAB method [56] and RNA contamination was removed using RNase. DNA concentration was checked via a Thermo Scientific NanoDrop and DNA integrity was assessed on a 0.8% agarose gel. Genotyping-bysequencing (GBS) was done following the published protocols [57]. After constructing GBS libraries as described by Alipour et al. [58], sequencing reads were trimmed to 64 bp and grouped into sequence tags, and SNP markers were called after alignment, which permits mismatches up to 3 bp. Markers were called in TASSEL software using the UNEAK pipeline. For avoiding false positive SNPs arising from sequencing errors, SNPs were filtered out if they had a missing rate > 10%, a MAF < 1%, and heterozygosity > 10%. The remaining missing was imputed using LD KNNi in TASSEL [58]. In the SNP calling pipeline, the wheat W7984 genome assembly was regarded as the reference genome [59].

Population structure and kinship matrix
The putative number of subpopulations (K) was determined by STRU CTU RE v2.2 using 10,000 burn-in iterations, followed by 10,000 proper MCMC sample steps for K-values ranging from K = 1 to K = 10 [60]. The best-fitting K value was determined using the ΔK method [61]. The matrix of population structure (Q) was calculated for the entire sample collection using a principal component analysis (PCA) implemented with the package Tidyverse in R. The kinship matrix (K) was obtained using the package GAPIT in R [62]. For cluster analysis, the elements of the kinship matrix were regarded as similarities and the outputs were visualized using UPGMA in GAPIT [63]. A neighbor-joining tree was constructed based on a pairwise distance matrix [63] and visualized by Archaeopteryx to determine the relationship between landraces and cultivars.

GWAS analysis
GWAS was carried out to detect marker-trait associations (MTAs) using the package mrMLM in R [21]. We considered −log 10 (P-value) ≥ 3.0 (P ≤ 0.001) as the significance threshold based on the previous reports [58,59]. All SNPs which met the above cut-off value were identified as significant MTAs. The GWAS results were visualized using Manhattan plots by the GAPIT package [64]. In the Manhattan plot, the x-axis and y-axis represent the chromosomal positions of SNPs and the −log 10 (P-value) is derived from the F-test, respectively. Q-Q plots were also obtained to further assess the results obtained from the Manhattan plots [23].

Candidate gene identification
To detect candidate genes affecting salinity tolerance during the seeding stage, regions surrounding traitsassociated SNPs were blasted against the rice and wheat genomes in the Ensemble genome database using the BLASTn. The IWGSC RefSeq v2.0 and IRGSP 1.0 were selected as genome references for wheat and rice, respectively [59,65]. After alignment, genes exhibiting the highest blast score and identity percentage were selected for gene ontology analyses.

Genomic prediction (GP)
The genomic prediction was performed using three different models: Bayesian ridge regression (BRR) [66], ridge regression-best linear unbiased prediction (RR-BLUP) [67], and genomic best linear unbiased prediction (GBLUP) [68]. All GP analyses were performed using the iPat software [69]. For three subpopulations, 10, 20, and 30% of genotypes were randomly assigned to a validation set with the remaining individuals used as the training set. For all of the GP procedures, the whole prediction process was repeated 100 times for each method. The accuracy of GP was presented as Pearson's correlation (r) between BLUPs and GEBVs over the training as well as validation sets.